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Abstract: Drought has pronounced and immediate impacts on agricultural production, especially in 
semi-arid and arid rainfed agricultural regions. Quantification of drought and its impact on crop yield is 
essential to agricultural water resource management and food security. We investigated drought and its 
impact on winter wheat (Triticum aestivum L.) yield in the Chinese Loess Plateau from 2001 to 2015. 
Specifically, we performed a varimax rotated principal component analysis on drought severity index (DSI) 
separately for four winter wheat growth periods: pre-sowing growth period (PG), early growth period 
(EG), middle growth period (MG), and late growth period (LG), resulting in three major subregional DSI 
dynamics for each growth period. The county-level projections of these major dynamics were then used to 
evaluate the growth period-specific impacts of DSI on winter wheat yields by using multiple linear 
regression analysis. Our results showed that the growth period-specific subregions had different major 
DSI dynamics. During PG, the northwestern area exhibited a rapid wetting trend, while small areas in the 
south showed a slight drying trend. The remaining subregions fluctuated between dryness and wetness. 
During EG, the northeastern and western areas exhibited a mild wetting trend. The remaining subregions 
did not display clear wetting or drying trends. During MG, the eastern and southwestern areas showed 
slight drying and wetting trends, respectively. The subregions scattered in the north and south had a 
significant wetting trend. During LG, large areas in the east and west exhibited wetting trends, whereas 
small parts in south-central area had a slight drying trend. Most counties in the north showed significant 
and slight wetting trends during PG, EG, and LG, whereas a few southwestern counties exhibited 
significant drying trends during PG and MG. Our analysis identified close and positive relationships 
between yields and DSI during LG, and revealed that almost all of the counties were vulnerable to 
drought. Similar but less strong relationships existed for MG, in which northeastern and eastern counties 
were more drought-vulnerable than other counties. In contrast, a few drought-sensitive counties were 
mainly located in the southwestern and eastern areas during PG, and in the northeastern corner of the 
study region during EG. Overall, our study dissociated growth period-specific and spatial location-specific 
impacts of drought on winter wheat yield, and might contribute to a better understanding of monitoring 
and early warning of yield loss. 
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1 Introduction 


Global climate change is one of the most complex challenges that human is facing in the 21* 
century. As recently indicated by the sixth assessment report of Intergovernmental Panel on 
Climate Change (IPCC, 2021), global warming as caused by human activities has led to more 
frequent and extreme climate events, such as heat waves and drought (particularly drought caused 
by increasing evapotranspiration). Drought has a pronounced and immediate impact on agricultural 
production, and consequently affects food security and social stability (Lesk et al., 2016; Madadgar 
et al., 2017; Kontgis et al., 2019). One of the key factors for better understanding and mitigating 
drought disasters is to assess agricultural drought and its impact on crop production. 

Agricultural drought can be defined as crop water deficit that results from two simultaneous 
processes: insufficient soil water supply caused by high temperatures or reduced precipitation, 
and continuous crop water consumption via transpiration. Consequently, the regular physiological 
processes of crops are not well-supported, thereby limiting crop growth or even leading to crop 
failure (Zhou et al., 2017; Dai et al., 2020). Researchers assess agricultural drought by various 
drought indices based on ground observations, remote sensing, or an integration of these two 
methods. Ground observation-based indices, such as the Palmer Drought Severity Index (PDSI; 
Palmer, 1965) and the Standard Precipitation Index (SPI; McKee et al., 1993) are limited by the 
sparse and heterogenetic spatial distribution of observation stations in agroecological regions, 
leading to insufficient support for timely drought detection, monitoring, and decision-making at a 
large spatial scale (Anderson et al., 2011; Liu et al., 2020). Remote sensing-based indices, such as 
the Vegetation Condition Index (VCI, Kogan, 1995) and the Temperature Vegetation Drought 
Index (TVDI; Sandholt et al., 2002) overcome some limitations of ground station-based indices, 
and potentially provide cost-effective and spatially continuous dynamic large-scale drought 
monitoring (Xu et al., 2018; Jiao et al., 2021). Evolving methodological developments have led to 
drought indices that integrate data from ground observations and remote sensing, such as the 
Vegetation Drought Response Index (VegDRI; Brown et al., 2008) and the Drought Severity 
Index (DSI; Mu et al., 2013). DSI integrates drought diagnostic information from the Normalized 
Difference Vegetation Index (NDVI) and the ratio of evapotranspiration (ET) to potential ET 
(PET), leading to flexible drought monitoring at different spatial and temporal scales. NDVI 
potentially links meteorological changes to vegetation responses as observed through changes in 
vegetation greenness with land-atmosphere water, carbon, and energy fluxes that are associated 
with climatic feedbacks (Atkinson et al., 2011). ET, as a key component of the terrestrial water 
and energy cycles, is an important constraint on water availability, and directly and effectively 
describes the moisture status of ecosystems (Yang et al., 2014; Srivastava et al., 2017). The ratio 
of ET to PET is commonly used as an indicator of terrestrial water availability and associated 
wetness or drought (Huang et al., 2020; Han et al., 2021). DSI was employed in this study for the 
assessment of drought because it permits more accurate characterization of drought indices while 
comprehensively considering these key factors at global and regional scales (Mu et al., 2013; 
Zhang and Yamaguchi, 2014; Zhang et al., 2016; Yu et al., 2019; Huang et al., 2020). 

Researchers can employ numerical models and statistical models to assess the impact of 
drought on agricultural production. Numerical models are process-based and can test the impact 
of water stress on crop yield using experimental trials that leverage a priori knowledge on crop 
physiology and reproduction, agronomy, and soil science (Lobell and Burke, 2010; Boonwichai et 
al., 2018; Lecerf et al., 2018; Bahri et al., 2019). Such models are not available for every region, 
and depend on site-specific calibration that requires a large amount of region-specific data such as 
information on cultivars, management, and soil conditions (Jizumi et al., 2009; Chapagain et al., 
2022). Alternatively, statistical models such as regression are flexible at both temporal and spatial 
scales (Mu et al., 2013). Despite their deficiencies in leveraging information about crop 
physiology and multicollinearity (Lobell and Asseng, 2017), they do provide several advantages 
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over numerical models. First, they require only the historical records of crop yield and 
corresponding drought metrics as previously mentioned (Mu et al., 2015); second, they are less 
dependent on field calibration data (Lobell and Burke, 2010); third, the uncertainties of the 
Statistical models can be transparently assessed (Lobell and Burke, 2010). Considering these 
advantages, we employed a statistical model based on multiple linear regression to assess the 
impact of drought on crop yield. 

Previous studies have found that China is a drought-prone country, and frequently suffers from 
drought-related crop failures (Wang et al., 2014; Geng et al., 2016), especially in the semi-arid 
and arid agricultural regions because croplands in such regions are mainly subject to natural 
precipitation (Hu et al., 2014; Qin et al., 2014; Bo et al., 2015). The Loess Plateau (LP) is a 
typical rainfed agronomic region in China, and is dominated by continental semi-arid and arid 
conditions. The agricultural production of this area is especially dependent on weather conditions, 
and is vulnerable to drought (Yu et al., 2019). The dominant crop in LP is winter wheat (Triticum 
aestivum L.). The growing season occurs over late-autumn, winter, and early-spring, and is 
characterized by insufficient precipitation (He et al., 2015; Jin et al., 2018). As a result, the yield 
of winter wheat in this region is typically relatively low and unstable (Jin et al., 2018). 
Consequently, this region is particularly suitable for investigating the impact of drought on 
agricultural production. 

Recent studies have revealed complex trends regarding the evolution and spatiotemporal 
characteristics of drought. Therefore, it may be necessary to assess drought at the subregional 
scale (Portela et al., 2017; Guo et al., 2018; Zhou et al., 2020). However, most previous studies 
have overlooked the evaluation of drought in subregions with similar drought co-variability when 
evaluating the impact of drought on winter wheat yield. In addition, water stress during different 
growth periods, such as the pre-sowing period before seeding and other growth periods after 
seeding, may influence crop phenology differently. Consequently, the relationship between 
drought and final yield should be evaluated separately for different crop growth periods. However, 
previous studies in LP either did not consider this effect or only considered growth periods after 
seeding (Zhang et al., 2016; Jin et al., 2018; Huang et al., 2020; Wan et al., 2021). Additionally, 
previous remote sensing-based studies regarding the impact of drought on crop yield have mainly 
focused on an entire crop production area (Lobell and Burke, 2010; Liu et al., 2018; Wang et al., 
2018a), with limited specificity and practical guidance for a specific type of crop (for example, 
winter wheat). 

Therefore, to address these issues, the general objective of this study was to assess wheat 
growth period-specific drought and its impact on yield in the Loess Plateau, with the 
consideration of the complex spatiotemporal characteristics of drought at the subregional scale. 
Specifically, we sought to extract reliable growth period-specific spatiotemporal dynamics of DSI. 
We also intended to dissociate growth period-specific and spatial location-specific impacts of 
drought on winter wheat yield. To this end, we first extensively assessed the subregional temporal 
dynamics of DSI by applying principal component analysis (PCA) with varimax rotation 
separately for different growth periods, and then examined the significant trends of DSI for each 
county according to the subregions and their corresponding major temporal characteristics in each 
growth period. We then evaluated the relationships between the major DSI dynamics during 
different growth periods with winter wheat yield individually for each county by using multiple 
linear regression analysis. We expect that our results will clarify the quantitative relationship 
between agricultural drought and winter wheat yield, and provide methodological considerations 
that can be applied to other regions for the evaluation of drought and its impact on crop yield over 
large-scale areas. 


2 Materials and methods 


2.1 Study area 
The Loess Plateau (LP; 33°43’—41°16'N and 100°54’—114°33’E) is located near the geographical 
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center of China, and covers about 64x10* km?, consisting of 334 county-level administrative 
regions that are scattered across seven provinces (Fig. 1). The region is characterized by a gradual 
transition from semi-arid conditions in the west to semi-humid conditions in the east, and 
experiences a typical continental monsoonal climate. The annual mean air temperature ranges 
from 3.6°C to 14.3°C, and the annual mean precipitation varies from 200 mm in the northwest to 
800 mm in the southeast (Fu et al., 2017). Precipitation fluctuates across seasons, and about 60% 
to 80% of the annual precipitation occurs from June to September in high-intensity rainstorms. A 
previous study indicated that from 2000 to 2012 annual ET had increased by 3.4 mm. However, 
during that same period, annual precipitation had only increased by 2.0 mm, suggesting a possible 
pressure on water availability (Jin et al., 2017). Such unbalanced seasonal precipitation results in 
the area suffering frequently from both soil erosion and drought (Jiang and Bai, 2020; Wu et al., 
2021). 

The LP is characterized by high solar radiation, relatively low precipitation, high 
evapotranspiration, and severe water deficit, and is mainly fed by rain (Li et al., 2010). 
Consequently, crop growth and yield are impacted by frequent drought occurrences. One of the 
major land cover types in the LP is cropland (Yu et al., 2019). We chose the main food crop in the 
LP, winter wheat, to evaluate the effects of drought on yield. Accordingly, the winter wheat 
producing area, about 2.51x10° hm?, was retrieved from a previous study (Jin et al., 2018) as our 
study area. 
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Fig. 1 Overview of the winter wheat producing counties in the Chinese Loess Plateau 


2.2 Data preparation 


2.2.1 Winter wheat yield 

Winter wheat production and sown areas from 2001 to 2015 were obtained from the China 
Statistical Yearbook and Provincial Statistical Yearbook for each county of the study area. We 
calculated the wheat yield of each county according to the corresponding production and sown 
area. 

2.2.2 Normalized difference vegetation index (NDVI) 

The NDVI dataset of the study area (1-km spatial resolution and 16-d temporal resolution) was 
acquired from the National Aeronautics and Space Administration (NASA) Earth Observation 
System (http://modis.gsfc.nasa.gov) from 2001 to 2015. We corrected the NDVI time series using 
a Savitzk-Golay filter to suppress the residual noise from cloud contamination and atmospheric 
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variability (Savitzky and Golay, 1964; Chen et al., 2004). Then the corrected data were 
interpolated to 1-d resolution using Lagrange interpolation (Bognar et al., 2017; Tian et al., 2019). 


2.2.3 ET 


We estimated actual ET and potential ET (PET) of winter wheat cropland from 2001 to 2015 
according to our previously published Vegetation Interfaces (VI)-based study (Wang et al., 2020) 
that involved the input of cropland cover maps, vegetation forcing data, and meteorological 
forcing data. A winter wheat cropland map was retrieved from the winter wheat producing area of 
a previous study (Jin et al., 2018). The vegetation forcing data included satellite-observed leaf 
area index (LAI), albedo, and solar zenith angle data. LAI and albedo were derived from datasets 
from the Global Land Surface Satellite (http://glass-product.bnu.edu.cn/) with an interval of 8 d 
and a spatial resolution of 0.05 degree. Solar zenith angle was acquired from the daily global 
0.05-degree resolution surface reflectance datasets from the Advanced Very High-Resolution 
Radiometer Long Term Data Record (https://www.ncdc.noaa.gov/cdr/terrestrial/avhrr-surface- 
refectance). The meteorological data included daily precipitation, air temperature, air pressure, 
relative humidity, sunshine duration, and wind speed. These data were obtained from the China 
Meteorological Administration (http://www.cma.gov.cn/) for 65 meteorological stations within 
and near the LP, and were interpolated to raster maps that covered the entire study area by using 
the gradient inverse distance square method that considers the effects of elevation (Nalder and 
Wein, 1998). All of the input data were first resampled to 1 km/?, then interpolated to 1-d 
resolution using Lagrange interpolation (Bognár et al., 2017), and finally input to the VI-based 
model to estimate four water vapor fluxes of winter wheat cropland on a daily basis, i.e., 
vegetation transpiration, evaporation from the soil surface, canopy interception, and PET. The 
actual daily ET was estimated as the sum of the first three of these water vapor fluxes, with 
relatively good agreement (Pearson correlation coefficient=0.69, NRMSE (normalized root mean 
squared error)=0.18, MAE (mean absolute error)=0.56 mm/d). More detailed methodological 
procedures for the VI-based model can be found in previous studies (Liang et al., 2020; Wang et 
al., 2020). 


2.3 Analysis of climatic yield anomalies 


Winter wheat yield is usually influenced by weather parameters and human activities such as 
improved crop varieties and advances in agricultural technology (Wang et al., 2018a). To 
eliminate the interference of human activities, we removed the yield trend from the wheat yields 
using the Hodrick-Prescott (HP) filter (Wang et al., 2014). This detrended yield, i.e., climatic 
yield, reflects the effect of climatic factors on the yields. We then standardized the climatic yield 
anomalies according to the following formula (Potopova et al., 2015): 


_ CY -Hoy 


Ocy 


AY (1) 


where AY is the climatic yield anomaly; CY is the climatic yield (kg/hm?); jy is the interannual 


mean of climatic yield (kg/hm7); and oy is the standard deviation of the climatic yield. 


2.4 Analysis of DSI 


DSI is an estimation of drought derived from NDVI and the ratio of ET to PET in which NDVI 
potentially links meteorological changes and vegetation responses (Atkinson et al., 2011), and the 
ratio of ET to PET is an index related to the water, carbon, and energy cycles of the land surface 
(Mu et al., 2013). DSI is capable of monitoring agricultural drought severity and predicting the 
effect of drought on crop yields (Zhang et al., 2016). This index was assessed from 2001 to 2015 
for the study area via four steps. The first step was to obtain the Z-score of the ratio of ET to PET. 
The second step was to obtain the Z-score of NDVI. The third step was to calculate a total Z-score 
represented by the sum of the Z-scores in the first and second steps. The last step was to acquire 
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DSI as the Z-score of the total Z-score. The calculation of a Z-score is a standardization process 
described as: 


Z, aA Mi (2) 


where x is the input data, which can be either the ratio of ET to PET, NDVI, or the sum of the 
Z-scores of the above two indices; 4» is the long-term mean of the input data; o, is the standard 
deviation of the input data; and Zy is the transformed Z-score of the input data. The further 
quantification of drought levels is frequently expressed as equal interval classification (Table 1) 
(Mu et al., 2013). 


Table 1 Categories for wet and dry conditions of the drought severity index (DSI) 


Category DSI Category DSI 
Extremely wet 21.50 Incipient drought —0.30 —-0.59 
Very wet 1.20-1.49 Mid drought —0.60 — -0.89 
Moderately wet 0.90-1.19 Moderate drought —0.90 — -1.19 
Slightly wet 0.60-0.89 Severe drought —1.20 — -1.49 
Incipient wet spell 0.30-0.59 Extreme drought <-1.50 
Near normal 0.29- —0.29 


To assess DSI across the winter wheat growing season at the county level, we divided the 
season into four growth periods according to previous studies: pre-sowing growth period (PG), 
early growth period (EG), middle growth period (MG), and late growth period (LG). PG was 
defined as the three-month period prior to sowing date. EG was the period from the sowing date 
to the green-up date; green-up dates were estimated from NDVI using a polynomial-cumulative 
NDVI algorithm (Liu et al., 2017; Wang et al., 2017). MG was regarded as the period from the 
green-up date to anthesis (generally two months later than the green-up date) (He et al., 2015). 
LG was the period from grain-filling to maturation (generally one month after anthesis) (He et al., 
2015). We assessed the spatiotemporal patterns of DSI across these growth periods. Specifically, 
we performed the following procedures for each growth period. We first applied a PCA to the 
yearly DSI time-series at the county level, resulting in a set of linearly uncorrelated principal 
components (PCs) and their mixing coefficients (Raziei et al., 2009). To get the main 
characteristics of DSI, we then selected a number of PCs in each growth period for further 
analysis such that at least 70% of the total variance of the original DSI data was expressed 
(Huang et al., 2015; Gupta and Jain, 2018). These selected PCs and their mixing coefficients were 
rotated using the Varimax method in order to stably express localized spatial patterns of drought 
(Richman, 1986; Huang et al., 2015). Then, for each rotated PC (RPC), we normalized the dot 
product of the variance of the RPC and rotated mixing coefficients (Guo et al., 2018), resulting in 
a normalized loading matrix. Subsequently, we set 0.6 as the threshold for the normalized loading 
matrices according to previous studies (Portela et al., 2017; Wang et al., 2018b), leading to 
subregions that each had a specific temporal pattern of DSI characterized by the corresponding 
RPC. Finally, we performed the Mann-Kendall rank correlation test (MK-test) on each RPC to 
examine the temporal pattern of subregional DSI (Mitchell et al., 1966). Each test generated three 
types of statistical parameters, i.e., Z-score, UB and UF. Z-score indicated the significance of the 
trend, and UB and UF time series indicated abrupt changes via their intersections (Shi and Wang, 
2015). The trend of each RPC was estimated using Sen's method (Sen, 1968; Shi et al., 2018). 


2.5 Evaluation of drought stress on winter wheat yield 


A multiple linear regression model was employed to evaluate the effect of drought on the yield of 
winter wheat during different growth periods (Jung et al., 2017; Wan et al., 2021). Specifically, 
for each county in a growth period, we reconstructed its main DSI temporal pattern according to 
the subregional RPC and corresponding rotated mixing coefficient. These main DSI temporal 
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patterns for the four growth periods were then used as regressors to evaluate their effects on the 
yield of that county, described as: 

Yao = BX, + 2X2 + BX; + p4X4 +c, (3) 
where Yac is the climatic yield anomaly of the county; X; (=1, 2, 3, 4) is the main DSI temporal 
pattern of a growth period; £; (i=1, 2, 3, 4) is the regression coefficient corresponding to X;; and c 
is the residual of the model. 


3 Results 


3.1 PCA-based division of subregions 

To examine the subregional patterns of drought, we employed PCA with varimax rotation on DSI 
during different growth periods. We selected the first three PCs during each growth period for 
further analysis, and found cumulative variances of 83.51%, 71.42%, 70.59%, and 71.27% for PG, 
EG, MG and LG respectively (Fig. 2). The selected PCs were rotated using the varimax method 
to acquire RPCs (rotated PCs). The cumulative variances for RPCs were the same as for PCs, 
indicating that the total variances of PCs were evenly expressed by RPCs (Fig. 2). 


100 


PC3 mmm PC2 mm PCI 


83.51 83.51 3 RPC3 mm RPC2 mm RPCI 
80F | 9.69 
15.23 71.42 71.42 70.59 70.58 71.27 71.26 
8.41 10.28 7.64 
60 20.42 


40 


Percentage of variance (%) 


Growth period 


Fig. 2 Percentages of variances explained by the first three principal components (PC1—PC3) for DSI during 
different growth periods of winter wheat in the Chinese Loess Plateau. PG, pre-sowing growth period; EG early 
growth period; MG, middle growth period; LG late growth period; DSI, drought severity index; RPC, rotated 
principal component. 


The normalized loadings of RPCs were visualized to examine the subregions during different 
growth periods (Fig. 3). During PG, RPC1 identified a subregion that covered a large area in the 
northwestern part of the study area (Fig. 3a), whereas the subregions for RPC2 and RPC3 covered 
a relatively smaller area in the south and east, respectively (Fig. 3b and c). During EG, the 
subregion for RPC1 was located mainly in the west-central area (Fig. 3d). The subregion for 
RPC2 was located mostly in the northeastern and western area (Fig. 3e). RPC3 affected mainly 
southeastern area (Fig. 3f). During MG, the subregions for RPC1 and RPC3 were primarily 
located in the east and southwest, respectively (Fig. 3g and i), whereas RPC2 was scattered in the 
north and south (Fig. 3h). During LG, the subregions for RPC1 and RPC3 covered a large area in 
the east and west, respectively (Fig. 3j and 1), while the dynamics of RPC2 were observed in only 
a small part of the central southern area (Fig. 3k). 


3.2 Temporal dynamics of DSI 


To assess the temporal dynamics of DSI in each subregion, we evaluated RPCs with the MK-test 
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Fig. 3 Spatial patterns of DSI in subregions of the Chinese Loess Plateau, defined by the first three rotated 
principal components (RPC1—RPC3) during PG (a-c), EG (d—f), MG (g-i), and LG (j-k). PG pre-sowing growth 
period; EG early growth period; MG, middle growth period; LG, late growth period. Dark red areas indicate 
subregions that were identified by setting 0.6 as the threshold for the normalized loading matrices. DSI, drought 
severity index. 


(Fig. 4). During PG, RPC1 presented a significant rapid increasing trend (slope=1.45; P<0.01) 
with slight fluctuation. The trend was more stable after 2008, as indicated by UF (Fig. 4a). RPC 2 
showed a slightly decreasing trend, but the change was not significant (Fig. 4b). RPC3 was in an 
unstable state of alternating between dryness and wetness (Fig. 4c). During EG, RPC1 and RPC3 
did not exhibit any significant trends (Fig. 4d and f), whereas RPC2 showed a significantly 
increasing trend (slope=0.50; P<0.01) (Fig. 4e). During MG, RPC] and RPC3 each showed 
insignificant decreasing and increasing trends (Fig. 4g and i), while RPC2 increased significantly 
(slope=0.83; P<0.01) over time (Fig. 4h). During LG, RPC1 showed a significant increasing trend 
(slope=0.94; P<0.01) (Fig. 4j), whereas RPC2 slightly decreased with an insignificant trend, and 
presented the driest year in 2013 (Fig. 4k). RPC3 exhibited an abrupt change in 2012, with large 
interannual fluctuations before 2011, but showed a significant increasing trend after 2012 (Fig. 
4l). 


3.3 Spatial patterns of DSI trends 


We evaluated the spatial patterns of DSI trends during each growth period based on previous 
subregion divisions and RPC analysis. During PG, we observed significant increasing trends 
covering most of the northern counties, as well as significant decreasing trends covering a few 
counties in the southwest (Fig. 5a). During EG, a significant but less strong increasing trend was 
found in counties located at the northern border (Fig. 5b). During MG, a few counties presented 
significant trends in which increases were scattered in the northwest and northeast, while 
decreases were found in the south (Fig. 5c). During LG, significant increasing trends were found 
in nearly all of the counties (Fig. 5d). 


3.4 Effect of DSI on winter wheat yield 


Given the spatiotemporal variability of DSI during each growth period, we expected spatial 
variation of its effect on winter wheat yield. In this respect, we employed a multiple linear 
regression analysis of winter wheat yield with the major DSI temporal dynamics of the growth 
periods as regressors. The regression coefficients were visualized in spatial maps to evaluate the 
coordinated variations between winter wheat yields and DSI (Fig. 6). In general, during different 
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Fig. 4 Temporal dynamics and Mann-Kendal (M-K) statistics of DSI (drought severity index) in each subregion 
during PG (a-c), EG (d—-f), MG (g-i), and LG G-k) of winter wheat in the Chinese Loess Plateau. PG, pre-sowing 
growth period; EG, early growth period; MG, middle growth period; LG late growth period. The Mann-Kendal 
test generated three types of statistical parameters, i.e., Z-score, UB and UF. Z-score indicated the significance of 
trend, and UB and UF indicated abrupt changes via their intersections. Slopes are estimated trends of RPCs. **, 
|Z|>2.58, significant trend at 99% confidence level. 


growth periods, dryness/wetness conditions showed very different impacts on winter wheat yield. 
The largest regression coefficients (absolute median=1.11) were mainly observed during LG, 
suggesting that dryness/wetness conditions were more closely related to winter wheat yields. PG 
(absolute median=0.64) showed a slightly higher value compared with EG (absolute 
median=0.59), indicating a slightly stronger effect of dryness/wetness condition on yields during 
PG than during EG. MG showed the lowest mean value (0.44), suggesting that the yields were 
less related to the dryness/wetness condition during this period. Specifically, during PG, negative 
covariations were mainly clustered in the northwest, while positive values occurred mainly in the 
south-central and eastern areas (Fig. 6a). During EG, more negative covariations were observed 
during PG (Fig. 6b). During MG, most of the covariations were slightly positive (Fig. 6c), while 
during LG, covariations for almost all of the study area were relatively high and positive (Fig. 6d). 


4 Discussion 


4.1 Growth period-specific spatiotemporal dynamics of DSI 


The effect of drought on winter wheat yield is a dynamic process. Drought stress affects winter 
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Fig. 5 Trends of DSI (drought severity index) during PG (a), EG (b), MG (c), and LG (d) of winter wheat in 
each county of the Chinese Loess Plateau. Counties shaded by dots showed significant trends in their main DSI. PG, 
pre-sowing growth period; EG, early growth period; MG, middle growth period, and LG, late growth period. 
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Fig. 6 Spatial patterns of regression coefficients between DSI (drought severity index) and winter wheat yield 
anomalies during PG (a), EG (b), MG (c), and LG (d). PG, pre-sowing growth period; EG early growth period; 
MG middle growth period; LG late growth period. 


wheat yield differently during different growth periods (Farooq et al., 2014; Li et al., 2018). 
Consequently, the quantitative assessment of such effects requires accurate divisions of growth 
periods, and relies on the accurate identification of growth period-specific spatiotemporal dynamics of 
drought. Based on these considerations, we first defined winter wheat growth periods. These periods 
considered both the pre-sowing and post-sowing growth periods because sufficient pre-sowing 
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soil moisture is crucial to the effective sowing, emergence, and overwintering of winter wheat (Li 
et al., 2017). In addition, we dynamically estimated the green-up date according to NDVI data 
using a polynomial-cumulative algorithm (Liu et al., 2017; Wang et al., 2017), and subsequently 
defined post-sowing periods as EG (sowing date to green-up date), MG (green-up date to anthesis) 
and LG (grain-filling to maturation). Afterwards, by applying a PCA with varimax rotation on 
DSI, we effectively identified three distinct subregions for each growth period, each 
corresponding to specific main DSI temporal dynamics. These subregional DSI dynamics were 
projected to the county level, from which we observed wetting trends in most northern counties 
during PG, EG, and LG, and found drying trends in some southern counties during PG and MG. 
A number of previous studies have documented climate change as one of the important factors 
affecting the dryness/wetness conditions in the LP. Yao et al. (2015) found that the change of 
terrestrial surface dryness and wetness conditions during the crop growth period in the LP were 
primarily driven by temperature, precipitation, and PET. These dryness/wetness conditions in the 
south and north varied in opposite directions, and were caused by the westerly wind in the north 
and the subtropical high in the south (Yao et al., 2015). Liu et al. (2021) showed that the annual 
mean temperature and precipitation presented obvious increasing trends after 2000, and that 
temperature was mainly influenced by the Atlantic multidecadal oscillation, while precipitation 
was mainly affected by the El Nifio-Southern oscillation and the Arctic oscillation. In addition, 
the dryness/wetness conditions of the study area were also affected by human activities. Wang et 
al. (2020) showed that the increasing ET in the LP was mainly driven by the ongoing optimization 
of human land-use management that was closely related to socio-economic development. Studies 
in the LP also found that the land use/land cover changes caused by the Grain for Green project 
modulated both surface humidity (Zhang et al., 2018) and vegetation responses to drought (Ding 
et al., 2021). 


4.2 Effect of drought on winter wheat yield 


The growth period-specific subregional DSI dynamics presented different levels of effects on 
winter wheat climatic yields. Our results indicated that winter wheat climatic yields were more 
positively related to DSI during MG and LG than during PG and EG. Particularly during MG, the 
regression coefficients in a large portion of the counties in the study area turned from negative 
values in the preceding period to positive values, suggesting that effective impact of drought 
appeared, and gradually increased from the wheat green-up date to the anthesis date, especially in 
the northeastern and north-central parts of the study area. This impact can be explained by the fact 
that water stress before anthesis reduces the number of grains per spike (Ji et al., 2010; Hlaváčová 
et al., 2018). Afterwards, during LG, the coefficients in most of the counties turned more positive, 
and were the highest across the growth periods, indicating an increased water requirement from 
grain-filling to maturation as compared with previous period. Indeed, the grain-filling rate and 
duration, partially constrained by the sufficiency of the water supply, will affect grain size, and 
therefore impact yield (Farooq et al., 2014; Wu et al., 2018). During PG and EG, relatively 
smaller areas presented positive regression coefficients between climatic yield and DSI. These 
drought-sensitive areas were mainly located at the southwestern and eastern parts of the study 
area during PG, and in northeastern part of the study area during EG. As for the remaining two 
growth periods, DSI was either not noticeably related to climatic yields or negatively co-variated 
with climatic yields, suggesting only small effects of drought on wheat climatic yields. This can be 
caused by multiple reasons. For example, the overwintering and green-up stages during EG are 
characterized by relatively low vegetation coverage, during which solar radiation is the dominant 
factor affecting winter wheat growth and development, while soil moisture plays a less important 
role (Huang et al., 2020). In addition, other agricultural experiments have shown that slight 
drought during EG increases root growth, thereby allowing the root system to absorb moisture 
better during subsequent growth stages (Wang et al., 2021). Conversely, waterlogging during this 
period may not only negatively impact yield, but may also increase the vulnerability of the crop to 
subsequent water deficit (Dickin and Wright, 2008). In summary, during PG, droughts in the 
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southwestern and eastern part of the study area may reduce winter wheat yields, while during EG, 
such drought-related yield reductions may appear only in several counties of the northeastern part 
of the study area. During MG, drought-vulnerable winter wheat production regions were mainly 
located at the northeastern and eastern part of the study area, while during LG such regions 
covered almost all of the study area. Our findings are generally consistent with a number of 
previous studies in the LP and other neighboring regions. For instance, a previous DSI-based 
study evaluated the impact of drought on winter wheat yield in Shaanxi, Henan, and Shanxi 
provinces, and found varied impacts on winter wheat yield across growing seasons, with the most 
significant impacts occurring during the heading and grain-filling stages (MG and LG) (Zhang et 
al., 2016). A study using multiple drought indices in northern China also found that winter wheat 
yield was sensitive to water stress during jointing and grain-filling stages (MG and LG). The 
authors explained this phenomenon in terms of the impact of water stress on the number of grains 
per spike, and suggested that crops can be protected from water stress during these two periods to 
maximize the probability of good harvests (Huang et al., 2020). Such remote sensing-based 
conclusions were also supported by a recent experimental study conducted at a site in the LP in 
which the authors observed modulations of winter wheat yield by drought stress during LG 
(Naseer et al., 2022). Taking all of these results together, we recommend that winter wheat should 
be protected from water stress mainly during MG and LG by rationally allocating the input of 
water and fertilizer in order to maximize the probability of a good harvest. 


4.3 Limitation and outlook 


Our current study had several limitations. First, DSI model relies on NDVI and ET/PET in order 
to improve drought monitoring. However, NDVI responds to drought with a time lag. A future 
study could further improve DSI by effectively considering this time lag in order to more reliably 
evaluate and monitor drought. Second, the dryness/wetness condition is normally considered as 
the dominant factor that impacts crop yields in semi-arid areas, but other factors, such as fertilizer, 
climatic anomalies, or economics also play roles in the agricultural production (Chen et al., 2015; 
Salehnia et al., 2020). The use of a single factor only permits the assessment of its impact on yield, 
rather than yield prediction. In a future study, we will attempt to evaluate and compare the 
contributions of other factors. Lastly, our study only considered linear relationships between DSI 
and yield, whereas complex nonlinear relationships may also exist (Lobell and Burke, 2010). A 
future analysis could also introduce nonlinear models, and attempt to interpret the potential 
nonlinear impacts. 


5 Conclusions 


The results of this study effectively and reliably identified growth period-specific subregional 
spatiotemporal patterns of drought in the winter wheat cropping area of the LP, China. Growth 
period-specific co-modulations between drought and winter wheat yield were found. Yield was 
most sensitive to drought during LG for almost the entire study area. Similar but less strong 
relationships existed during MG mainly over the northeastern and eastern counties. There were 
few drought-sensitive counties during PG and EG. These findings contribute to a better 
understanding of spatiotemporal relationships between drought and winter wheat yield in the LP. 
Overall, our study demonstrates a novel methodological approach for dissociating growth 
period-specific and spatial location-specific impacts of drought on a specific type of crop, which 
may contribute to a more effective monitor and forecast of drought-induced yield loss, and may 
further provide scientific support for better guiding of a sustainable development of regional 
agricultural management. 
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